
%%%%%%%%%%%%%%%%%%%
%%% Simulation %%%
%%%%%%%%%%%%%%%%%%%

clear all
tic

for i=[3 5 7] % three different seeds

%% Load Calibration Outcome
load calibration_outcome.mat

%% Setup Model Parameters

    params.alpha=0.5;
    params.eta=6;
    params.r=0.02;
    params.pi=0.001;

    params.lambdaE=xx_cal(1);
    params.lambdaS=xx_cal(2);
    params.lambdaI=xx_cal(3);
    params.lambdaC=xx_cal(4);
    params.f=xx_cal(5);
    params.xiE=xx_cal(6);
    params.xiI=xx_cal(7);
    params.lambdae=params.lambdaE;
    params.xie=xx_cal(8);

%% Setup Equlibrium Variables

    eqm.xI=xx_eqm(1);
    eqm.xE=xx_eqm(2);
    eqm.xe=xx_eqm(3);
    eqm.gq=xx_eqm(4);

%% Setup seed
    seed=i*100

%% Setup Algorithm Parameters

  algS.delt      = 0.005;                % for converting continuous time to discrete time: implies that a year has 1/delt subperiod
  algS.nFirm     = 2^18;                 % the number of starting firms
  algS.maxFirm   = floor(1.1*algS.nFirm);% max number of firms
  algS.nBurn     = 100;                  % Number of years for burning
  algS.nReal     = 10;                   % Number of years for calculating statistics 

  algS.NNExt     = 50;                   % Max. number of product lines
  
%% Run Simulation

[nProd,age,cohort,firmID,quality,ageProd,cohortProd,ProdID,numBS,numCN,check] = simulation_function(params,eqm,seed,algS);

% pre-allocation
combined=zeros(size(quality,1)*size(quality,2),7,size(quality,3)-1);

% combining all product/firm information
for q=1:algS.nReal*4
    ready_firmID=repelem(firmID(:,q+1),1,algS.NNExt);
    ready_age=repelem(floor(age(:,q+1)/(1/algS.delt)*4),1,algS.NNExt);
    ready_cohort=repelem(floor(cohort(:,q+1)/(1/algS.delt)*4),1,algS.NNExt);
    ready_ProdID=ProdID(:,:,q+1);
    ready_quality=quality(:,:,q+1);
    ready_ageProd=floor(ageProd(:,:,q+1)/(1/algS.delt)*4); % in quarters (keeping integer)
    ready_cohortProd=floor(cohortProd(:,:,q+1)/(1/algS.delt)*4); % in quarters (keeping integer)
    combined(:,:,q)=[ready_firmID(:) ready_age(:) ready_cohort(:) ready_ProdID(:) ready_quality(:) ready_ageProd(:) ready_cohortProd(:)];
end
    combined_ready=[combined(:,:,1) ones(size(combined,1),1)];
    data_for_Stata=combined_ready((combined_ready(:,4)>0),:);
for q=2:algS.nReal*4
    combined_ready=[combined(:,:,q) q*ones(size(combined,1),1)];
    data_for_Stata=[data_for_Stata;combined_ready((combined_ready(:,4)>0),:)];
end
    data_for_Stata=num2cell(data_for_Stata);
    filename = sprintf('../STATA/INPUT/data_for_Stata_seed%d00.csv', i);
    dlmwrite(filename,data_for_Stata,'precision',8)
    clear all
end

toc
